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THE SUPERMASSIVE BLACK HOLE AND DARK MATTER HALO OF NGC 4649 (M60) 
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ABSTRACT 

We apply the axisymmetric orbit superposition modeling to estimate the mass of the supermassive 
black hole and dark matter halo profile of NGC 4649. We have included data sets from the Hubble 
Space Telescope, stellar, and globular cluster observations. Our modeling gives M, = (4.5 ± 1.0) x 
1O 9 M and M/Ly- obs = 8.7 ± 1.0 (or M/L y = 8.0 ± 0.9 after foreground Galactic extinction is 
corrected). We confirm the presence of a dark matter halo, but the stellar mass dominates inside 
the effective radius. The parameters of the dark halo are less constrained due to the sparse globular 
cluster data at large radii. We find that in NGC 4649 the dynamical mass profile from our modeling is 
consistently larger than that derived from the X-ray data over most of the radial range by roughly 60% 
to 80%. It implies that either some forms of non-thermal pressure need to be included, the assumed 
hydrostatic equilibrium may not be a good approximation in the X-ray modelings of NGC 4649, or 
our assumptions used in the dynamical models are biased. Our new M, is about two times larger 
than the previous published value; the earlier model did not adequately sample the orbits required to 
match the large tangential anisotropy in the galaxy center. If we assume that there is no dark matter, 
the results on the black hole mass and M/Ly obs do not change significantly, which we attribute to 
the inclusion of HST spectra, the sparse globular cluster kinematics, and a diffuse dark matter halo. 
Without the HST data, the significance of the black hole detection is greatly reduced. 
Subject headings: black hole physics - galaxies: general — galaxies: nuclei — galaxies: kinematics and 
dynamics — stellar dynamics 



1. INTRODUCTION 

Most nearby galaxies harbor supermassive black holes 
at their centers. Correlations bet ween black hole 
(BH) mass and host galaxy p roperties (iMagorrian et al.l 
19981: IGebhardt et all 120001: iFerrarese fc MerrittJ |2000F 
Haring fc Rixll2004D have been used extensively in theo- 
retical models in order to understand growt h of the black 
hole and galaxy (e.g., iHopkins et alJl2"008| ). The latest 
work from Hopkins et al. suggest that a main role of a 
black hole is to halt star formation in the galaxy when the 
black hole is la rge enough , thereby causing the dichotomy 
in colors (e.g., Bell 2008). While there is still an active 
debate as to the relative role of AGN feedback versus 
star formation feedback, there is a consensus that physi- 
cal mechanisms for black hole growth are very important 
for understanding mass growth in galaxies. A concern is 
that the black hole correlations may have significant sys- 
tematic biases, both from kinematics with poor spatial 
resolution and models that do not adequately include th e 
full mass profile (see discussion in IGultekin et all[2009h . 

Dynamical modeling of galaxies using orbit super- 
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models do not limit the form of the allowed velocity 
anisotropics. Thus, the stellar orbital structure resulting 
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from the dynamical modeling provides a unique window 
on the mass growth process in the massive systems, as 
long as the model assumptions are valid. 

A particular sy st emat ic bias is shown in 
IGebhardt fc Thomas) (|2009D where they find that 
the black hole mass can be underestimated in the most 
massive galaxies if the dark halo is not included. They 
find a degeneracy between the dark halo and black hole 
mass, since without the dark halo the stellar mass-to- 
light ratio is overestimated which subsequently decreases 
the required contribution of the bl ack hole to constrain 
the ce ntral kinematics. In M87, [Gebhardt fc Thomasl 
(2009) find that the black hole mass goes from 2.5 x 10 9 
to 6.4 x 10 9 Mq by simply running models including 
a dark halo. Furthermore, the uncertainties do not 
overlap, implying a large systematic effect. In M87, 
however, this degeneracy is strong since the black hole 
is not well resolved by the kinematic data. For galaxies 
with well spatially-resolved kinematics, we do not expect 
the degeneracy to be as significant. 

In addition to studying the black hole mass, there is 
a strong need to study the shapes of dark matter pro- 
files. There is still little consensus fo r those measured in 
indiv idual galaxies (e.g. . PNe from iRomanowskv et al.l 
12001 stel lar light from iForestell fc Gebhardt! 120081 X- 
rays from iGastaldello est alJl2007t iChurazov et alJ 12008b 
i Humphrev et al.l I2008L 120091. and globu lar clusters from 
iBridges et afl !2006: H wang et al J 120081 ) . The impressive 
work using gravitational le nsing to measure the avera ge 
dark matter profiles (e.g., Mandclb aum et all I2006allb1 ) 
has been able to reach out to nearly 1 Mpc. However, 
these results need to be compared to measurements based 
on individual galaxies. It is important to understand the 
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galaxy-to-galaxy scatter in the profiles and whether there 
are environmental effects. 

This paper is part of an extensive campaign to mea- 
sure both the black hole mass and the dark matter pro- 
file simultaneously, to examine the possible bias of the 
dark matter profile on the black hole mass estimate for 
galaxies with various profiles, and to compare with the 
gravitational potentials derived with other independent 
techniques such as the X-rays and weak lensing studies. 
Initially, we focus on the more massive galaxies, but it 
will be impo rtant to extend future analy sis over a large 
mass range. IGebhardt fc Thomas! (|2009t ) report results 
for M87, the most massive nearby elliptical, and in this 
paper we focus on NGC 4649. Both galaxies are giant 
ellipticals with a central surface brightness "core" . 

In this paper we present the axisymmetric orbit su- 
perposition models for NGC 4649 (M60), combining 
data from the Hubble Space Telescope (HST), stellar, 
and globular cluster observations. NGC 4649 is a gi- 
ant elliptical with low surface brightness located in a 
subclump to the e ast of the main Virgo concentration 
(|Forbes et al.ll2004T ). NGC 4649 has been studied exten- 
sively for its total mass profile in recent X-ray model- 
ings (e.g.. Humphrey et al.ll2006l : iGastaldello et all 120071 : 
Humphrey et al.ll2008D . and globular cluster studies (e.g., 
Bridges et al.1 12006'). The goals of our study are to place 
NGC 4649's black hole mass estimate on a more solid 
footing, to infer the properties of its dark matter halo, 
and more importantly to offer an independent cross- 
check using different dynamical tracers to the previous 
studies on NGC 4649. 

We assume a distance to NGC 4649 of 15.7 Mpc. At 
this distance, 1" corresponds to 76 pc. 

2. DATA 
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Fig. 1. — The V-band volume luminosity density profile for stars. 
This comes from the deprojection o f the surface brightness profile 
presented in Kormcndv et al. (2009). 

The input data of NGC 4649 (surface brightness pro- 
file and kinematics) have all been previously published, 



and we rely on that data for the models presented here. 
W e use the stellar surfac e brightness profile compiled 
by iKormendv et al.1 (|2009f) from a v ariety of sources, in - 
cluding HST i maging as present e d in | Lauer et al.1 {2005). 
All data from IKormendv et al.1 (|2009( ) have been trans- 
formed to y-band, and extend to about 700". We de- 
project the surf a ce bri ghtness using Abel's formula as in 
IGebhardt et al.l ((199 6.) . The volume luminosity density 
we derived is plotted in Figure [U The effective radius 
R e of NGC 4649 is about 90" by fitting a de Vaucouleurs 
profile to its surface brightness profiles. 

T he stellar kin e matics are the same as used 
in IGebhardt et ail (|2003l), and are presented in 
iPinknev et all ((20031 ). These data include long-slit ob- 
servations on HST STIS and th ree position angles fr om 
ground-based observations (see Pink nev et al.1 120031 for 
details). There are no new stellar kinematic data pre- 
sented in this paper. The HST spectra extend to 0.8", 
and the ground-based stellar data extend to 70". The 
dynamical modeling code uses the line-of-sight velocity 
distribution (LOSVD) , which are the same as used in 
IGebhardt etaD IpOOl . 

In order to extend to yet further radii, we also include 
data from globular cluster (GC) velociti es. The GC kine- 
matics come from IHwang et al.l (120081) . whic h combine 
kinem atic data from iBridges et aL I (|2006l) and Hie et al.l 
((2001) . There are 121 clusters with velocities over a ra- 
dial range from 32-533". However, we only use the ve- 
locities that are beyond 200", since inside that radius the 
clusters with velocities are spread over too many model 
bins to add any significant information compared to the 
stellar kinematics. A further complication with using 
GCs in the central parts is that their number density 
profile is not well known there. Outside of 200", there 
are 61 clusters with velo cities, and the ave rage radius 
of these clusters is 319" ((Hwang et al.l[2008l ). Rotation 
in the GC system at radii be yond 200" is reporte d in 
IHwang et all (|2008fk however, IBridges et all ((20061 ) ar- 
gue for little to no rotation at radii around 200". We 
therefore explore models where we include the rotation 
directly (by fitting a velocity profile on the major axis 
at 319" with a rot ation amplitude an d velocity disper- 
sion as reported bv lHwang~et al. 20081) and where we use 
the second moment of the GC velocities only (by fitting 
one velocity profile, centered on zero velocity, to the GC 
velocities in the full annulus). The results are nearly 
identical. In this paper, we only present results from us- 
ing the second moment. Our analysis of the Hwang et al. 
data implies a second moment of 267(±25) km/s , whic h 
is consistent with the analysis of Bridg eset al.l {2006). 
Since the modeling code uses LOSVDs directly, we must 
transform the moments into LOSVDs, which we do using 
Monte Carlo re alizations of the value s and uncertainties 
(as discussed in IGebhardt et al.1l2003| ). 

The dynamical models require the number density pro- 
file of the GCs (since we use clusters as a tracer pop- 
ulation). IHwang et al.l ((20081) present number density 
profiles of the GCs. We use their equation 5 for the com- 
bined sample as the deprojected number density profile. 
Comparing this profile to the stellar luminosity density 
shows that beyond 150", the two profiles have the same 
slope. Inside that radius, the GC number density profile 
flattens more compared to the stars. Since we do not use 
any clusters inside of 200", we cadopt the same number 
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density profile for the clusters as we do for the stars. 

3. DYNAMICAL MODELS 

The orbit supe rposition models we u se are based on the 
idea invent ed bvlSchwarzschildl (11979D. an d are described 
in det ail in IGebhardt et al l (120001 l200l ; iThomas et all 
(|200l [200l ; lSiopis et al.l (|2009f TThe models have been 
tested extensively to recover both the black hole mass 
and the dark matter halo without significant bias for ax- 
isymmetric systems. 

We first determine the luminosity density from the sur- 
face brightness profile as described in Sj2] Although we 
have constructed models using different inclinations, we 
generally assume that the galaxy is edge-on and adopt 
the minor-to-major axial ratio of 0.9. Next we com- 
pute the gravitation field from the 3-D light distribu- 
tion assuming that the mass consists of a black hole 
with mass M m , and stars with a constant M/Ly obs ratio 
independent of position, and a dark matter halo (now 
parametrized with only logarithmic and NFW profiles) . 

Next our orbit-based model is implemented in two 
steps: constructing the orbit library and fitting the full 
LOSVD (as opposed to using only the second veloc- 
ity moments) with the orbit library. We measure the 
LOSVD in radius and angle relative to the symmetry 
axis of the galaxy, and in velocity. We use 24 radial, 5 
angular, and 15 velocity bins. The radial range of the dy- 
namical models extend to 1000". We calculate the galaxy 
potential and force on a grid that is 5 times finer than the 
grid used in the data comparison. We carefully choose 
the initial positions of the orbits to generate a sufficientl y 
dense sampling of E, L z , and I3 ([Thomas et al.l [20041) . 
This ensures a reliable representation of the full phase 
space. Generally we have at least 30, 000 orbits when 
both the black hole and dark matter halo are included. 

We then choose the orbital weights so the superposi- 
tion of orbits match the light distribution and LOSVDs 
of the galaxy as well as possible. The orbital weights are 
derived from the maximization of the objective function 
S = S — ax 2 where S is an approximation to the usual 
Boltzmann entropy and x 2 is the sum of squared resid- 
uals to the data. The smoothing parameter a controls 
the influence of the entropy S on the orbital weights. We 
cannot specify the optimal value of a in advance. In prac- 
tice we start with a very small a and gradually increase 
it until the improvement in x 2 is less than a given per- 
centage in a single iteration. We have verified that our 
results are not sensitive to the choice of the smoothing 
parameter a. 

We ran nearly 16,000 different models to estimate the 
black hole mass and the dark matter halo parameters. 
Each model has a distinct set of orbit library, and takes 
about 1.5 hours on the supercomputer lonestar of the 
Texas Advanced Computing Center. 

4. RESULTS 

4.1. Models including a dark matter halo 

In our canonical set of models, the mass distribution p 
of NGC 4649 consists of a central black hole, stellar mass 
density, and a dark matter halo: 

p = M,5(r) + Tv + pbm 



where T = M/L v obs is the mass-to- light ratio of the 
stars (a position- independent constant), and v is the 
observed T^-band stellar luminosity density. M/L v obs 
needs to be corrected with the foreground Galactic ex- 
tinction in order to compare with that derived from stel- 
lar population models. A p opular foreground G alactic 
extinction value A v = 0.088 ([Schlegel et al-l fl^S! could 
be adopted; it gives the extinction-corrected M/L v = 
M/L v>obs /1.084. 

We describe the dark matter halo with a logarithmic 
profile, whose density is given as 



Pdm 



V? 3r 2 



AttG (r 2 +r 2 ) 2 



and potential as 



$ = -y c 2 ln(l + -) 



IGebhardt k Thomas! ([2009D found that realistic NFW 
models that fit the globular cluster kinematics are not 
centrally concentrated enough to dominate the mass in 
the inner regions; the enclosed mass profile they get with 
NFW models are very similar to those obtained with 
logarithmic halos. Thus, whether we parameterize with 
an NFW or a logarithmic halo will have little influence 
on the black hole mass estimate, particularly true for 
NGC 4649 where HST spectra exist. We focus on a log- 
arithmic halo in this study. 

The kinematics include all data sets from the HST, 
stellar, and globular cluster observations We use 

the x 2 distribution from all possible models to determine 
the four best-fit parameters (M/L y obs , M„ V c , r c ), and 
the associated uncertainties. We start with a uniform but 
sparse grid in the 4-D parameter space, then sample the 
smallest x 2 region with a finer grid. As we will see from 
the figures below, the x 2 minimum and the contours of 
NGC 4649 are quite regular, so our sampling procedure 
should be adequate to cover the parameter space near the 
best-fit values. The uncertainties in the parameters are 
determined from the change in the marginalized x 2 as we 
vary one of the variables; A% 2 = 1 above the minimum 
represents the 68% confidence band or la uncertainty. 
We use the middle of the 68% confidence range as the 
best-fit value and half of that range as the uncertainty. 

Figure [2] presents % 2 as a function of M/L y obs , the 
black hole mass M., the halo scale velocity V c , and the 
halo core radius r c , including all values for the other three 
parameters. Each point represents a possible model, but 
we show only models near the x 2 minimum to highlight 
the la uncertainty. The solid line along the bottom ridge 
represents the marginalized x 2 values which we use to de- 
termine the best fit and uncertainties. Given the number 
of LOSVDs, the velocity binning used in the modeling, 
an d that the LOSVD b ins are correlated (as discussed 
in IGebhardt et all 120031 ). the reduced x 2 is around 0.4, 
which is typical with orbit-based models. 

From Figure[2j we nn d the black hole mass M. = (4.5± 
1.0) x 10 9 M Q and the stellar M/L^ obs = 8.7 ± 1.0 (or 
M/Ly = 8.0 ±0.9 after foreground Galactic extinction is 
corrected with Ay — 0.088). The scale velocity V c and 
the core radius r c of the logarithmic halo are not well 
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Fig. 2. — x 2 versus M/L v obs , black hole mass M,, the halo scale velocity V c , and the halo core radius r c . Each point represents a 
possible model for the logarithmic halo (only models with A\ 2 < 25 over the minimum are shown in order to highlight the la uncertainty). 
The solid line along the bottom ridge represents the marginalized \ 2 values which we use to determine the best fit and uncertainties. The 
horizontal line marks the 68% confidence limit (Ax 2 = 1.0). 



constrained. As we will discuss later, halo parameters 
are degenerate with each other, and their exact values 
have little effect on the shape of the dark matter profile 
in the region of interest. 

In order to find the possible correlations and degen- 
eracies in the four parameters, Figure |3] plots M/h v obs , 
M,, V c and r c against each other. The small grey points 
represent the locations of all the models. The large black 
and red points highlight those models that arc within the 
68% (Ax 2 = 1.0) and the 95% (Ax 2 = 4.0) confidence 
bands, respectively, after marginalizing over the other 
possible parameters. The strong correlation between V c 
and r c is apparent. This is because the dynamical mod- 
els only fit for the enclosed mass and we do not have 
enough kinematics at large radii; we cannot distinguish 
dark halos with a large V c and a large r c from the one 
with a small V c and an accordingly fine-tuned small r c . 
By the same token, the degeneracy between M/L v obs 
and M, (apparent in the the 95% confidence band) is ex- 
pected; the decreasing/increasing contribution from the 
stars can be made up by having a larger/smaller black 
hole mass, to certain extent. For NGC 4649, the degener- 



acy between M/L v obs and M. is relatively weak because 
the inclusion of HST data makes t he uncertainties in M, 
much smaller, unlike in M87's case ([Gebhardt fe Thomasl 
2009). HST spectra provide good spatial sampling inside 
of the region where the black hole dominates over stars. 
We do not find obvious degeneracies among other param- 
eters. 

The radial profile of the enclosed mass in NGC 4649 
is shown in Figure 2] The black lines represent the 
models that are within the 68% confidence band of the 
best fit (large black points in Figure [3]) . The red line 
stands for the stellar contribution alone using the best- 
fit M/L v obs of 8.7. The stellar component clearly dom- 
inates the mass profile from 10" to 100". At the ef- 
fective radius R e (~ 90"), the best-fit stellar mass ac- 
counts for about 75% of the total enclosed mass. Within 
10" the 4.5 x 10 9 M Q black hole is dominant, while out- 
side 100" the dark matter halo prevails. Compared with 
NGC 4649 the best-fit dark halo of M 87 is much more 
dominant ( Ge bhardt fc Thomasl [20091 ) . This also shows 
that ellipticals with comparable luminosities can have 
quite different distributions of dark matter relative to 
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FlG. 3. — Plots of the V-band M/Ly obs , black hole mass Mm (in units of 10 Mg), the scale velocity V c of the logarithmic halo in km/s, 
and the core radius r c of the logarithmic halo in kpc against each other to show possible degeneracies. The small grey points represent the 
locations of all the models. The large black points are the locations of those models that are within the 68% confidence limit (Ax 2 = 1.0) 
after marginalizing over the other possible parameters. The red points are those that are within the 95% confidence band. The strong 
correlation between V c and r c is apparent. 
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Figure [5] presents the integrated mass-to-light ratio 
M/Ly obs as a function of radius. The black dotted lines 

are the total M/Ly obs of the models that are within the 
68% confidence band of the best-fit potentials (black lines 
in Figure [4J . The black solid lines are the same models 
but excluding the contribution from the black hole. Con- 
sistent with Figure _H Figure [5] shows that the increase 
of the total M/Ly obs at small radius is dictated by the 
black hole, and the increase at large radius is due to the 
dark matter halo. In the intermediate radial range, from 
10" to 100", the M/Ly obs is a constant that stellar con- 
tribution dominates. 

The shaded band in Figure [5] shows the wide range 
of possible M/L y obs from stellar population models de- 
rived as below. We first adopt the age (« 12 Gyr) 
and metallic i ty \Z /H] « 0.25 for NGC 4649 from 
iTrager et al.l ([2000). With the age and metallicity, we 
then interpolate for the stellar M/Ly from tablefl pro- 

2 http:/ /www. icg. port. ac.uk/~maraston/SSPn/ml/ml_SSP.tab 



duced by the evolutio nary population synthesis models 
(IMarastodll9^[2T)05l ). If the Salpeter IMF (0.1 to 100 
Mq) is assumed, then we get the interpolated stellar 
M/Ly = 8.06 as the upper bound of the shaded area. If 
the Kroupa IMF (0.1 to 100 Mq) is used instead, then we 
obtain the interpolated M/Ly = 5.13 as the lower bound. 
In order to compare with our uncorrected M/Ly obs in 
Figure [5j we have to multiply the M/Ly from stellar 
population models by a factor of 1.084 (Ay = 0.088) to 
mimic the foreground Galactic extinction. From Figured] 
we see that our best fit M/Ly obs = 8.7 ± 1.0 agrees best 
with the stellar population modeling result assuming the 
Salpeter IMF, but it is also consistent with the range of 
results assuming different IMFs. 

4.1.1. Comparison with the X-ray results 

In Figure |4] we also include the the mass profile in- 
ferred from X-ray m odelings (Gastaldcllo et al.l 120071 : 
iHumphrev et al.ll2008l ). shown as the green circles. They 
are derived from the traditional approach, i.e., from the 
smoothed parametric fits to the temperature and density 
of the X-ray emitting gas assuming hydrostatic equilib- 
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Fig. 4. — The mass profile for NGC 4649. The black lines repre- 
sent the models that are within the 68% confidence band of the 
best fit (as in Figure [3}- The green line is our representation 
of the X-ray derived enclosed mass profile (the green circles, as 
in IGastaldello et al.l [20071 ; IHumphrev et al.l 120081 ) . The parame- 
ters for the green line are M. = 3.5 x 10 9 M , M/h v obs = 5.0, 

V c = 410km/s, and r c = lOkpc. The red line is the average con- 
tribution from the stars, where we integrate the light profile in 
Figure [T] then times 8.7 (the best-fit M/L v obs ). Clearly the total 
mass profile obtained from our modeling is consistently larger than 
that from X-rays over most of the radial range. 
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Fig. 5. — The integrated Nl/h v obs as a function of radius. The 

black dotted lines represent the models that are within the 68% 
confidence band of the best-fit potentials (black lines in Figure [1). 
The black solid lines are the same models but excluding the con- 
tribution from the black hole. The green solid and dotted are our 
best representation of the X-ray potential without and with a black 
hole. The shaded horizontal band is the range of possible M/L^ obs 

(after applying the foreground Galact ic extinction of A y = 0.088) 
from stellar population models (Maraston 1998, 2005), depending 
on what form of IMF is assumed. 

rium. The green line in Figure [4] is our best represen- 
tation of the X-ray derived mass, assuming a logarith- 
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Fig. 6. — x 2 distribution when a logarithmic halo with V c = 
410km/s and r c = lOkpc, which best matches the results of X-ray 
modeling (the green line in Figure l4t . is imposed in the modeling. 

(a) : x 2 as a function of black hole mass M m , including all values in 
M/L^ obs . Each point represents a model, however, not all mod- 
els that we ran are shown in order to highlight the la uncertainty. 
The solid line along the bottom ridge represents the marginalized 
X 2 values which we use to determine the best fit and uncertainties. 
The horizontal line marks the 68% confidence limit (Ax 2 = 1.0). 

(b) : as in (a) but for x 2 as a function of V-band M/L v obs , includ- 
ing all values in M. . (c) : 2-D plot of x 2 as a function of black hole 
mass and M/L v obs . The points represent models that we ran. As 

in [Gebhardt et al. (2002), the contours were determined by a 2-D 
smoothing spline interpolated from the these models and represent 
Ax 2 of 1.0, 2.71, 4.0 and 6.63 (corresponding to 68%, 90%, 95%, 
and 99% for 1 degree of freedom) . 

mic dark halo. The parameters for the green line are 
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Fig. 7. — The radial profile of Ax , which is the difference in 
the radius-cumulated x between the best-fit model (M. = 4.5 X 
10 9 M Q , M/h v obs = 9.0, V c = 550km/s, and r c = 15kpc) and 
the model that best matches the X-ray derived potential (M. = 



3.5 x 10 9 M , M/L 



V, obs 



5.0, V c = 410km/s, and r c = lOkpc). 



M. = 3.5 x 10 9 M , M/L y obs = 5.0, V c = 410km/s, 
and r c = lOkpc. In Figure [5] the green solid and dotted 
lines are the M/L v obs of our best representation of the 
X-ray potential without and with a black hole, respec- 
tively. 

The total mass profile obtained from our modeling is 
consistently larger than that from X-rays over most of 
the radial range by a factor of about 1.7. Our stellar 
mass alone is larger than the X-ray derived mass, so the 
difference cannot be due to the different dark matter ha- 
los derived from the two methods. Our best-fit black 
hole mass is consistent within the la error, yet about 
35% larger than of the M. = 3.35 x 1O 9 M derived in 
Humphrey et all ([2008). Detailed discussions on the dif- 
ference from the X-ray mass profile are in i j5.ll 

A concern from this work could be that the dark halo 
parameters are poorly constrained by the globular cluster 
kinematics, leading to a biased result on the mass pro- 
file and potentially on the black hole mass. We therefore 
include the X-ray data in two different ways to compare 
with our dynamical models. First, we use the X-rays to 
determine the large radii mass profile, allowing our dy- 
namical models to determine the stellar M/L and black 
hole mass. Secon d, we adopt the full X- ray mass pro- 
file as measured by Humphr ey et al.l (|2008[ ) and compare 
with our best-fit model. 

First, to use the large radii mass profile as determined 
from the X-rays, we run a grid of models assuming the X- 
ray dark halo potential (V c = 410km/s, and r c = lOkpc), 
and fit for black hole mass and stellar M/L v obs . In other 
words, we fix only the outer part of the green line in Fig- 
ure IU and allow the inner mass profile to seek the best fit 
to the data. The results are shown in the distribution of 
X 2 as function of M, and M/L y obs in Figure |U We still 
find very similar M. and M/Ly obs as in the canonical 
set of models ( ^4.ip ; both the best-fit values and the la 
uncertainties of M, and M/L v obs in Figure [5] are very 
similar to those found in Figure [5] and Figure [3J 

Second, we run the mass model with the X-ray derived 
parameters (A/. = 3.5 x 10 9 M , M/L v obs = 5.0, V c = 
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Fig. 8. — As in Figure \E\ but for the x 2 distribution when no 
dark halo is included in the modeling. Large radii kinematics are 
excluded in the models. The absolute values of x 2 in this figure 
are not comparable to those in Figure [2] and [6] as this figure fits to 
less kinematic data (75 data points from the one LOSVD for the 
globular kinematics and 4 LOSVDs for stellar kinematics beyond 
40" are excluded). 

410km/s, and r c = lOkpc), and find the best fit to the 
kinematic data. The \ 2 of this model is about 98, which 
is a much worse fit than our best-fit model (% 2 w 43.5). 
So the stellar and globular cluster kinematics are indeed 
inconsistent with the X-ray derived mass profile. Fig- 
ure [7] demonstrates how Ax 2 changes with radius, where 
A% 2 is the difference in the radius-cumulated x 2 between 
the best-fit model in Figure [2] (M. = 4.5 x 10 9 M Q , 
M/L v obs = 9.0, V c = 550km/s, and r c = 15kpc) and 
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the model with the X-ray derived parameters (the green 
line in Figure H M. = 3.5 x 10 9 M & , M/L^ obs = 5.0, 
V c = 410km/s, and r c = lOkpc). The roughly linear in- 
crease of Ax 2 with radius indicates that data from all 
radial range contribute equally to the difference in the 
total x 2 - This is another reason why the assumed outer 
potential from X-rays does not affect the best-fit M. and 
M/L y!0bs . 

4.2. Models without a dark matter halo 

iGebhardt fc Thomas! ()2009f ) found that including a 
proper dark matter halo in the dynamical modeling of 
M87 can give a factor of two larger black hole mass than 
excluding one. Here we also do a similar test on how 
the black hole mass and M/L v obs change if we exclude 

a dark halo in the modeling of NGC 4649 (i.e., the mass 
distribution is composed of the black hole and the stars 
only, p = M,5(r) + Tv). Since we know beforehand 
that the dark halo dominates at very large radii, we in- 
tentionally exclude some of the large radii kinematics in 
order to not be unduly influenced by data that extend 
well into the dark halo regime. The excluded regions are 
the one LOSVD for the globular cluster kinematics and 
the four LOSVDs for the stellar kinematics beyond 40". 
We test the effect of including all data, without a dark 
halo, and find that the resultant stellar M/L v obs (a con- 
stant independent of radius) is biased to a large value 
of 11.5, and the minimal \ 2 is larger by ~ 25. Thus, 
the large radii kinematics are excluded in the no-halo 
modeling. Figure [8] presents the x 2 distribution when 
a dark halo is excluded and the large radii kinematics 
are excluded. From the marginalized x 2 as a function 
of M, and M/L v obs in Figured! we find the black hole 

mass M. = (4.3 ± 0.7) x 10 9 M Q and the F-band stellar 
M/L v obs = 9.0 ± 1.0. The uncertainties with the no 
dark halo models are slightly lower than when including 
a dark halo, and this difference could be due to noise in 
the x 2 contours in Figure El The best- fit black hole mass 
is quite differe nt from the previously published ones for 
NGC 4649 in IGebhardt et all (|2003h . which are cruder 
models without a dark halo either; we will discuss this 
point in detail in $5] Also these values differ from those 
when a dark matter halo is allowed (§ I4.1[) by less than 
5%. So in NGC 4649 the inclusion of a proper dark 
halo does not make as big a difference to the best-fit M, 
and M/Ly obs as in the case of M87. The existence of 
HST data in NGC 4649 certainly helps to alleviate the 
bias. Also the globular cluster kinematics in NGC 4649 
are sparse (only one radial bin), so they do not have as 
much leverage in the total x 2 as m M87. Furthermore, 
the dark matter halo of NGC 4649 is not as massive and 
concentrated as that of M87, so the black hole mass es- 
timate of NGC 4649 is not biased as much by neglecting 
the dark matter halo. 

We further experimented recomputing for the best-fit 
model using the ground-based spectral data only, exclud- 
ing the HST spectra. The inclusion of the HST spectra 
in the modeling undoubtedly improves the significance 
of the black hole substantially; the % 2 difference between 
models with the best-fit M, and with zero M, is around 
30. The best-fit black hole mass is almost the same, how- 
ever, the uncertainty in M, becomes much larger when 



the HST spectra are excluded. This experiment empha- 
sizes the need for the HST spectra to improve the error 
in the determination of M t . 

4.3. Velocity Dispersion Tensor 
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Fig. 9. — Shape of the velocity dispersion tensor as a function 
of radius for the best-fit model. The black line is along the major 
axis, and the red line is near the minor axis. 



We also examine the internal orbital structures by 
studying the shape of the velocity dispersion tensor. In 
our work the tangential dispersion is defined as at — 

\J a 4> a e > wnere a 4> includes contributions from both 
random and ordered motions, i.e., it is the second mo- 
ment of the azimuthal velocity relative to the systemic 
velocity instead of to the mean rotation velocity. An non- 
rotating isotropic model gives o r /ot — 1. Figure [9] shows 
the internal dispersion ratio <J r /o~t along the major and 
minor axes for our best fit model (M, = 4.5 x 10 9 Mq, 
M/L y obs = 9.0, V c = 550km/s, and r c = 15kpc). Over 
most of the radial range, the model is very isotropic. It 
is strongly tangentially biase d near the center , and more 
radially biased at large radii. IGebhardt et all (|2003D dis- 
cussed the possibility that the amount of the tangen- 
tial anisotropy at large radii may be overestimated when 
the dark matter halo is excluded. The current result 
shows the orbits are clearly more radially dominated at 
r > 100" when the dark matter starts to dominate. 

We also want to inspect if the distribution function of 
NGC 4649 truly depends on three integrals of motion or 
only two: energy E and the z-componcnt of the angular 
momentum L z . We know that in a two-integral model, 
/ = f(E, L z ) = f(vl/2 + vl/2 + v 2 J2 + $(i?, z), R V<t> ), 
o~r must be equal to a z . One way to test whether a two- 
integral model is sufficient is by comparing ay and oq on 
the equatorial plane where an = a r 1^=0 and a z = erg | ,£=o 
(Figure fTTJj) . The radial variation of a r /a$ in Figure ITTJ1 
demonstrates that the best-fit model is inconsistent with 
having only two integrals of motion, despite that the 
radially-averaged <r r /ae is close to one. Thus, our three- 
integral model is needed to best fit the data. We also 
show the (70 /erg ratio in the bottom panel of Figure [TDJ 
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Fig. 10. — Ratios of the velocity dispersions cr r /&g (top) and 
u^/crg (bottom) as a function of radius near the equatorial plane 
for the best-fit model. 

We can see that near the black hole the tangential dis- 
persion is dominated by the azimuthal <f> motion, whereas 
the 9 and <fr dispersions are similar at most of radii. Such 
detailed information on the shape of velocity dispersion 
tensor could be potential constraints on the formation 
processes of the supermassive black holes. 

5. DISCUSSIONS 
5.1. Estimates of the mass profile 



From 34.1.11 we see that the total mass profile from 
our modeling is consistently larger than that from X- 
rays over most of the radial range by a factor of about 
1.7. Also the slopes of the two mass profiles inside the 
effective radius are quite similar. If our orbit-based mod- 
eling results are confirmed in other independent studies, 
then it may suggest that the X-ray derived mass is sys- 
tematically lower. There are several possible reasons why 
the X-ray modeling could underestimate the true mass 
profile. 

X-ray modelings generally assume that the X- 
ray emitting hot ga s is in hydrostatic equilibrium. 
Dic hl fc Statlerl ([2007) argue that the X-ray gas may not 
be in good hydrostatic equilibrium because they found 



no correlation between X-ray and optical ellipticities in 
the inner region where stellar mass dominates over dark 
matter and a shallow correlation should be expected. If 
the hot gas is not in hydrostatic equilibrium, the inflow 
of gas can mak e the X-ray derived mass smaller than 
the t rue value (jPellegrini fc CiottH [2006; Jo hnson et al.l 
I2009T) . 

Even if the hydrostatic equilibrium holds reason- 
ably well for the X-ray emitting gas, the existence 
of non-thermal pressure components, but omitted in 
the X-ray modeling, can still make the X-ray derived 
mass estimate less than the true value. For exam- 
ple, possible forms of non-thermal pressure support in- 
clud e magnetic field, cosmic rays , and microturbulence 
(e.g. IBrighenti fc Mathewi Il997t iChurazov etaLl 120081: 
IJohnson et al.H2009f) . AlsoT since the hot gas in ellipticals 
comes mainly from the stellar mass loss, infall or merg- 
ers, the gas could carry significant amount of angular mo- 
mentum. The gas can settle into rotationally supported 
systems as it flows in and spins up. The spherically- 
determined X-ray mass may again be an underestimate 
since the ro tation is ignored. 

However. IBrighenti et all (|2009[ ) study NGC 4649 and 
conclude that turbulant motion in the X-ray gas is not 
strong enough to bias the derived potential. Our dynami- 
cal results are in disagreement with the X-ray derived po- 
tential, possibly suggesting that other non-thermal pres- 
sures are present. While it is still unclear what combi- 
nation of systematics contribute to the difference in the 
mass estimates of the X-ray and orbit-based methods, 
it is imperative to extend our orbit-based modeling to a 
larger sample of galaxies that have been analyzed with 
X-rays, and to examine if the difference between the two 
methods is general. 

5.2. Black Hole Mass 

The black hole mass reported here is 4.5(±1.0) x 
10 9 M(T). Using nearly the same input data, 
IGebhardt etaLl (|2003(1 report a black hole mass of 
2. 0(±0.5) x 10 9 Mp, (no te the reported black hole mass 
of IGebhardt et all 120031 has two changes applied to it: 
first, we decrease the reported mass by 7% for the differ- 
ent distances, and second, we apply a 9% increase in the 
mass due to a nu merical error in th e previous calcula- 
tion as reported in Siopis ct al. 2009). While this is only 
a 1.5-ct difference, since we are using the same data, we 
should get the same result. Thus, any difference must 
be in the modeling code changes, the input assumptions 
(e.g., dark halo mass) or a combination. We find that 
the difference is due to a more complete orbit sampling, 
discussed below, but we have also checked a variety of 
other possibilities. 

Comparing Figure [5] (no dark halo) with Figure [5] (in- 
cluding a dark halo) shows that there is little difference 
in the black hole mass when including a dark halo. This 
is expected since the black hole's effect on the kinematics 
is well resolved with the HST kinematics. Other effects, 
such as the radial range used in the modeling (in the 
present paper we must extend the galaxy model out to 
larger radii to include the dark halo), appear to have an 
insignificant effect as well. 

The main modifications used for this study is to in- 
crease the orbit library substantially and to sample the 
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Fig. 11.— The modeled and observed LOSVDs at R = 0.18" on 
the major axis. The y-axis re presents the relative projected light 
for this spatial bin. The data llPinknev et al. 2003) are shown by 
the solid dots and their 68% confidence bands (the velocity bins 
at the edge have a large uncertainty because the original LOSVD 
extractions from Pinkncy ct al. 2003 do not extend that far in ve- 
locity space). The solid line is the LOSVD from the best-fit model 
in this pa per, and the dashed lin e is the LOSVD from the best-fit 
model in Gebhardt ct al. (2003). The current model is a signif- 
icantly improved fit due to the inclusion of the highest angular 
momentum orbits (i.e., nearly circular orbits) that can match bet- 
ter the flat-topped LOSVD. 

phase space differently. For this study, we ran two sets 
of models, one with 30,000 orbits and one with 13,000 
orbits. The results are nearly identical in terms of A\ 2 
(which is what we use to determine the values and un- 
certainties), however the minimum % 2 is lower by about 
30. iGebhardt et~aT1 (|2003l ) use 7000 orbits, with a x 2 
that is larger by 80. Given that we do not see a differ- 
ence in the current models using different orbit numbers, 
we are confident that orbit number is not a concern for 
NGC 4649. 

We do, however, find that the change in orbit sampling 
is the the main factor causing difference in the black hole 
mass of NG C 4649. The d i fferenc e in orbit sampling is 
presented in iThomas et al.l (|2004h , where we now use a 
sampling based on the density in the meridional plane 
wher eas the previous version (used in Gebhard t et al.l 
I2003D relied on sampling along the zero-velocity curve. 
In particular, the new sampling covers more completely 
the phase space occupied by the highest angular momen- 
tum orbits for a given energy, like near-circular orbits and 
shell orbits. 

Core galaxies appear to have si gnificant tangential or - 
bital anisotropy in their centers ([Gebhardt et al.| [2003). 
Thus, if the orbital coverage poorly sampled important 
regions of phase space occupied by tangential orbits, 
this may bias the best-fit model. Indeed, the orbital 
anisotropy is strongly biased towards tangential motion 
in the central regions (Section 14.31). much more s o tha n 
in the original orbit-models of Geb hardt et al.1 ([2003). 
In Figure [11] we exemplify the comparison of LOSVDs 
from the old and new orbit sampling, and the observed 



LOSVD. In the figure the HST observed LOSVD at 
R = 0.18" on the major axis is quite flat-topped. The 
modeled LOSVD in the current best-fit model (solid 
curve) is abl e to fit the data v ery w ell, whereas the best- 
fit model in [Gebhardt ct al] ([20031 ) (dashed curve) has 
trouble matching the profile. The poor fit from the pre- 
vious model is due to the lack of nearly circular orbits. 

Increasing the amount of tangential anisotropy causes 
the projected dispersion to drop in the central regions; 
thus, in order to match the observed projected velocity 
dispersion, the black hole mass has to be increased to 
compensate. This concern needs to be tested on other 
ga laxies, specifi c ally c ore galaxies. We note that the work 
of ISiopis et all ([2009) compare black hole masses from 
analytic models using the same orbital sampling pre- 
sented here. They find no bias in the black hole masses, 
but the central density of their galaxy (NGC 4258) is 
cuspy instead of cored. Galaxies with such a large core 
as NGC 4649 and large tangential orbital bias should be 
further studied. 

With M, = 4.5 x 10 9 Af Q and the effective stellar ve- 
locity dispersion a e — 385 km/s, N GC 4649 's position is 
above the latest M, — a e relation ([Gultekin et al.l [20091 ) 
by about 0.33 dex, but still within the intrinsic scat- 
ter of the relation. Increasing the largest black hole 
masses by a factor of 2 — NGC 464 9 in this paper and 
M87 in IGebhardt k Thomasl ([2009D — will have impor- 
tant c onsequences for und erstanding the upper mass end 
(as in lLauer et al.l l2007| ) and comparison with masses 
as estimated for quasars. There has been a long stand- 
ing problem as to why some quasars have black hole 
mass estimates approaching 10 10 M Q , whereas none have 
ever been measured that high even though the volumes 
surveyed should be large to see these; a factor of 2 to 
3 increase appears to resolve this issue. Furthermore, 
an underestimate of the black hole mass seriously ef- 
fects the physical correlations, used extensively to quan- 
tify the role of bla ck holes on galaxy evolution (e.g., 
Hopk ins et al.l2008l) , the amount of deviation fr om a sim- 
ple M. — cr* power-law relation ([Wvithd 1 20061). and the 
number density of the lar gest black holes ([Lauer et al~l 
[20071 [Beniardi et alJl2007T ) . 

5.3. Main uncertainties 

The uncertainties presented in this paper are statistical 
uncertainties from the \ 2 analysis. In order to translate 
a given change in y 2 into a significance (for example, we 
report the 68% confidence band based on A\ 2 = 1), we 
make two major assumptions. 

First, we assume that we understand the uncertain- 
ties of the kinematics and any correlations between the 
measurements. Correlations in the spectral extractions 
of the LOSVDs may cause biases in the determination 
of the parameter uncertainties such as the black hole 
mass (discussed in Houghto n et al.ll2006[ ). For the anal- 
ysis used in this paper, there are likely correlations due 
to smoothing o f the L O SVD, as d i scuss ed originally in 
IGebhardt etaTI ([2000D . iGebhardtl (|2004l ) tries to quan- 
tify this effect by running Monte Carlo simulations, start- 
ing from the noise in the spectra (and then running dy- 
namical mo dels for each realiza tion). A most complete 
analysis is in lSiopis et al.l ([20091 ) who find a similar result 
(see their Fig. 17). The uncertainty as measured from 
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the Ax 2 analysis is similar to, but smaller than, the un- 
certainty as measured from the spread in the best-fitted 
values from the realizations. This work suggests that the 
X 2 analysis provides realistic uncertainties. Yet, we know 
there are correlations in the LOSVDs — these correlations 
exist in the non-parametric analysis that was used in this 
paper but also exist in the b asis f unction approach ad- 
vocated by iHoughton et al.l (|2006|) due to degeneracies 
with the continuum placement for the spectra. In this 
analysis of NGC 4649, we rely on A% 2 = 1 reflecting the 
68% confidence band. From the tests we describe above, 
we argue that the statistical analysis is understood, al- 
though expande d ana lysis along the line s presented in 
iGebhardtl (|2004D and iSiopis et al.l ((2009) is greatly de- 
sired. We feel the more immediate need is to control 
systematic effects, discussed next. 

Second, our reported uncertainties do not include ef- 
fects from systematic uncertainties. As is clear from the 
change in the black hole mass from the previous study, 
there are systematic uncertainties that need to be ex- 
plored. These uncertainties will influence the black hole 
mass, orbital structure, and dark halo profile. We discuss 
above the concern with the orbital sampling. We note, 
again, that the orbit sampling issue is likely only impor- 
tant for galaxies that have an extreme orbit al structure 
as seen here in NGC 4649. IGebhardti (l200l tested the 
effect of orbit number on a power-law galaxy and found 
no bias in the black hole mass estimate. Running ana- 
lytic models with a range of orbital distributions would 
be worthwhile, as well as studying additional core galax- 
ies. For core galaxies, another important systematic is 
that they may be better modeled as triaxial as opposed to 
axisymmetric as assumed here (e.g. , van den Bosch et all 
120081 ). Ivan den Bosch fc de Zeeuwl ((20091) show that the 
black hole mass in the core galaxy NGC 3379 increases if 
the model is allowed to be triaxial. This increase could 
be particular to the intrinsic orientation for NGC 3379 
but continued studies including triaxiality is warranted. 
Other model assumptions such as the assumed inclina- 
tion, variation of the stellar mass-to-light ratio with ra- 
dius (as derived from stellar population models), and 
the assumed shape of the dark halo profile, for example, 



could be additional sources of systematic uncertainties. 



6. CONCLUSIONS 

We model the dynamical structure of NGC 4649 us- 
ing the high resolution data sets from HST, stellar, and 
globular cluster observations. Our main new results are: 

1. Our modeling gives M. = 4.5 ± 1.0 x 10 9 M© and 
M/L v obs = 8.7±1.0. Our new M. of NGC 4649 is about 
a factor of 2 larger than the previous result. We find that 
the earlier model did not adequately sample the orbits 
required to match the large tangential anisotropy in the 
galaxy center.. 

2. We confirm the presence of a dark matter halo in 
NGC 4649, but the stellar mass dominates inside the 
effective radius. The parameters of the dark halo es- 
pecially the core radius are less constrained due to the 
sparse globular cluster data at large radii. 

3. Unlike in the case of M87, the black hole mass from 
the dynamical modeling is not biased as much by the 
inclusion of a dark matter halo, because high-resolution 
HST spectra are available for NGC 4649, the globular 
cluster kinematics are sparse, and the halo is not as dom- 
inant inside the effective radius R e as that of M87. 

4. We find that in NGC 4649 the dynamical mass pro- 
file from our modeling is consistently larger than that de- 
rived from the X-ray data over most of the radial range 
by about 70%. It implies that either some forms of non- 
thermal pressure need to be included, the assumed hy- 
drostatic equilibrium may not be a good approximation 
in the X-ray modeling, or the assumptions used in our 
dynamical modeling create a bias. 
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